3. Example Solow model

In this jupyter notebook we will specify, solve and analyse a simple Solow model in ModelFlow.

#Required packages
import pandas as pd

# Modelflow modules
from modelclass import model
   
#for publication 
latex=0
model.widescreen()

3.1. Specify the model

We start by defining the logic of the Solow model in the Business Logic Language.

fsolow = '''\
Income          = a  * Capital**alfa * Labor **(1-alfa) 
Consumption     = (1-saving_rate)  * Income 
Investment      = Income - Consumption 
diff(Capital)   = Investment-Depreciation_rate * Capital(-1)
diff(Labor)     = Labor_growth * Labor(-1) 
Capital_intensity = Capital/Labor 
'''

3.2. Create a model class instance

After defining the Business Logic Language and storing it in the variable ‘fsolow’, we create a class instance called msolow.

msolow = model.from_eq(fsolow,modelname='Solow model')

Above the equation for Capital and Labor on the left hand side of the = (equal to) consist of an expressing diff(Capital) and diff(Labor). The equations are not normalized.

To solve a model in modelflow all equations has to be normalized. Meaning that the left hand side only consist of variables not expressions. So the function model.from_eq will normalize the equations as the first step before the model can be solved.

In this case first diff(Capital) is transformed to \(\Delta capital = capital-capital(-1)\). Then the lagged variables is moved to the right side of the =. The same goes for diff(labor).

So the normalized business language of the model now looks like:

msolow.print_model
FRML <> INCOME          = A  * CAPITAL**ALFA * LABOR **(1-ALFA)  $
FRML <> CONSUMPTION     = (1-SAVING_RATE)  * INCOME  $
FRML <> INVESTMENT      = INCOME - CONSUMPTION  $
FRML <> CAPITAL=CAPITAL(-1)+(INVESTMENT-DEPRECIATION_RATE * CAPITAL(-1))$
FRML <> LABOR=LABOR(-1)+(LABOR_GROWTH * LABOR(-1))$
FRML <> CAPITAL_INTENSITY = CAPITAL/LABOR  $

3.3. Create some data

To show what Modelflow can do, we create a Pandas dataframe with input data. And print the first 5 out of 300 observations.

N = 300  
df = pd.DataFrame({'LABOR':[100]*N,
                   'CAPITAL':[100]*N, 
                   'ALFA':[0.5]*N, 
                   'A': [1]*N, 
                   'DEPRECIATION_RATE': [0.05]*N, 
                   'LABOR_GROWTH': [0.01]*N, 
                   'SAVING_RATE':[0.05]*N})
df.head(2) #this prints out the first 5 rows of the dataframe
LABOR CAPITAL ALFA A DEPRECIATION_RATE LABOR_GROWTH SAVING_RATE
0 100 100 0.5 1 0.05 0.01 0.05
1 100 100 0.5 1 0.05 0.01 0.05

3.4. Run the model

result = msolow(df,keep='Baseline') # The model is simulated for all years possible 
result.head(29)
LABOR CAPITAL ALFA A DEPRECIATION_RATE LABOR_GROWTH SAVING_RATE INVESTMENT CAPITAL_INTENSITY INCOME CONSUMPTION
0 100.000000 100.000000 0.5 1.0 0.05 0.01 0.05 0.000000 0.000000 0.000000 0.000000
1 101.000000 100.025580 0.5 1.0 0.05 0.01 0.05 5.025580 0.990352 100.511609 95.486029
2 102.010000 100.076226 0.5 1.0 0.05 0.01 0.05 5.051924 0.981043 101.038487 95.986562
3 103.030100 100.151443 0.5 1.0 0.05 0.01 0.05 5.079029 0.972060 101.580575 96.501546
4 104.060401 100.250762 0.5 1.0 0.05 0.01 0.05 5.106891 0.963390 102.137821 97.030930
5 105.101005 100.373733 0.5 1.0 0.05 0.01 0.05 5.135509 0.955022 102.710176 97.574667
6 106.152015 100.519926 0.5 1.0 0.05 0.01 0.05 5.164880 0.946943 103.297593 98.132713
7 107.213535 100.688931 0.5 1.0 0.05 0.01 0.05 5.195002 0.939144 103.900030 98.705029
8 108.285671 100.880357 0.5 1.0 0.05 0.01 0.05 5.225872 0.931613 104.517449 99.291576
9 109.368527 101.093830 0.5 1.0 0.05 0.01 0.05 5.257491 0.924341 105.149813 99.892323
10 110.462213 101.328993 0.5 1.0 0.05 0.01 0.05 5.289855 0.917318 105.797092 100.507238
11 111.566835 101.585506 0.5 1.0 0.05 0.01 0.05 5.322963 0.910535 106.459257 101.136294
12 112.682503 101.863045 0.5 1.0 0.05 0.01 0.05 5.356814 0.903983 107.136282 101.779468
13 113.809328 102.161300 0.5 1.0 0.05 0.01 0.05 5.391407 0.897653 107.828145 102.436738
14 114.947421 102.479976 0.5 1.0 0.05 0.01 0.05 5.426741 0.891538 108.534829 103.108087
15 116.096896 102.818793 0.5 1.0 0.05 0.01 0.05 5.462816 0.885629 109.256317 103.793501
16 117.257864 103.177483 0.5 1.0 0.05 0.01 0.05 5.499630 0.879920 109.992597 104.492967
17 118.430443 103.555792 0.5 1.0 0.05 0.01 0.05 5.537183 0.874402 110.743661 105.206478
18 119.614748 103.953478 0.5 1.0 0.05 0.01 0.05 5.575475 0.869069 111.509502 105.934027
19 120.810895 104.370310 0.5 1.0 0.05 0.01 0.05 5.614506 0.863915 112.290118 106.675612
20 122.019004 104.806070 0.5 1.0 0.05 0.01 0.05 5.654275 0.858932 113.085509 107.431233
21 123.239194 105.260550 0.5 1.0 0.05 0.01 0.05 5.694784 0.854116 113.895678 108.200894
22 124.471586 105.733554 0.5 1.0 0.05 0.01 0.05 5.736032 0.849459 114.720631 108.984599
23 125.716302 106.224895 0.5 1.0 0.05 0.01 0.05 5.778019 0.844957 115.560378 109.782359
24 126.973465 106.734397 0.5 1.0 0.05 0.01 0.05 5.820747 0.840604 116.414931 110.594185
25 128.243200 107.261893 0.5 1.0 0.05 0.01 0.05 5.864215 0.836394 117.284305 111.420090
26 129.525631 107.807224 0.5 1.0 0.05 0.01 0.05 5.908426 0.832323 118.168518 112.260093
27 130.820888 108.370242 0.5 1.0 0.05 0.01 0.05 5.953380 0.828386 119.067591 113.114212
28 132.129097 108.950808 0.5 1.0 0.05 0.01 0.05 5.999077 0.824578 119.981548 113.982470

3.5. Create a scenario and run again

dfscenario = df.upd('LABOR_GROWTH + 0.002')  # create a new dataframe, increase LABOR_GROWTH by 0.002
scenario   = msolow(dfscenario,keep='Higher labor growth ') # simulate the model 

3.6. Now the results are also embedded in msolow.

  • .basedf contains the first run of the model

  • .lastdf contains the last run of the model

Also in this case the keyword keep is used. This causes the results to be stored in a dictionary msolow.keep_solutions. This can be useful when comparing several scenarios.

3.7. Inspect results

3.7.1. Using the [ ] operator

We can select the variables of interest with wildcards. This will operate the results stored in basedf and .lastdf

3.7.1.1. Look at variables starting with a C

msolow['#ENDO']

3.7.1.2. Look at all endogenous variables

msolow['labor*'].dif.plot() 
../../../_images/Example Solow_23_0.png ../../../_images/Example Solow_23_1.png

3.7.2. Using the keept solutions

As mentioned above, because the keyword keep was used. The results are also stored in a dictionary. These data can also be used for charting.

The reason for placing the results in a dictionary is to enable comparison of many scenarios, not just the first and the last.

with msolow.set_smpl(1,30):
    msolow.keep_plot('income con*' ); 
../../../_images/Example Solow_25_0.png ../../../_images/Example Solow_25_1.png
msolow.modeldash('INCOME',jupyter=1)
apprun
Dash app running on http://127.0.0.1:5001/
Dash_graph(mmodel=<
Model name                              :          Solow model 
Model structure                         :         Simultaneous 
Number of variables                     :                   11 
Number of exogeneous  variables         :                    5 
Number of endogeneous variables         :                    6 
>, pre_var='INCOME', filter=0, up=1, down=0, time_att=False, attshow=False, all=False, dashport=5001, debug=False, jupyter=1, show_trigger=False, inline=False, lag=False, threshold=0.5, growthshow=False)

3.8. More advanced topics

3.8.1. The logical stucture

Now the model has been analyzed, and the structure can be displayed.

You will find more on the logical structure here

3.8.1.1. Model structure

msolow.drawmodel( title="Model Structure", png=latex,size=(15,15))
../../../_images/Example Solow_30_0.svg

3.8.1.2. Adjacency matrix

Another way to illustrate the dependency graph is an adjacency matrix.

msolow.plotadjacency();
../../../_images/Example Solow_32_0.png

The variables [‘INVESTMENT’, ‘CONSUMPTION’, ‘CAPITAL’, ‘INCOME’] in the red area are the core of the model and has to be solved as a system.

LABOR is the prolog and can be calculated before the core is solved. While CAPITAL_INTENSE is the epilog which can be calculated after the core is solved.

Many models comform to this pattern. And for solving purpose a model is divided into a prolog, core and an epilog. Even if the core is actually consistent of several strong components.

3.8.2. The python function used to solve the model

In order to solve the model Modelflow will generate a python function which implements the model. The user will hopeful newer have to relate to the generated python code. The point of modelflow is, that the user has to relate to the specification of the business logic, not the implementation in code

print(msolow.make_los_text)
def make_los(funks=[],errorfunk=None):
    import time
    import tqdm
    from numba import jit
    from modeluserfunk import jit, recode
    from modelBLfunk import array, classfunk, clognorm, exp, gamma, inspect, jit, lifetime_credit_loss, log, logit, logit_inverse, lognorm, matrix, mv_opt, mv_opt_prop, norm, normcdf, qgamma, sqrt, sum_excel, transpose
    def prolog0(values,outvalues,row,alfa=1.0):
        try :
            pass
            values[row,0]=values[row-1,0]+(values[row,5]*values[row-1,0])
        except :
            errorfunk(values,sys.exc_info()[2].tb_lineno,overhead=9,overeq=0)
            raise
        return 
    def prolog(values,outvalues,row,alfa=1.0):
        prolog0(values,outvalues,row,alfa=alfa)
        return  
    def core0(values,outvalues,row,alfa=1.0):
        try :
            pass
            values[row,9]=values[row,3]*values[row,1]**values[row,2]*values[row,0]**(1-values[row,2])
            values[row,10]=(1-values[row,6])*values[row,9]
            values[row,7]=values[row,9]-values[row,10]
            values[row,1]=values[row-1,1]+(values[row,7]-values[row,4]*values[row-1,1])
        except :
            errorfunk(values,sys.exc_info()[2].tb_lineno,overhead=20,overeq=1)
            raise
        return 
    def core(values,outvalues,row,alfa=1.0):
        core0(values,outvalues,row,alfa=alfa)
        return  
    def epilog0(values,outvalues,row,alfa=1.0):
        try :
            pass
            values[row,8]=values[row,1]/values[row,0]
        except :
            errorfunk(values,sys.exc_info()[2].tb_lineno,overhead=34,overeq=5)
            raise
        return 
    def epilog(values,outvalues,row,alfa=1.0):
        epilog0(values,outvalues,row,alfa=alfa)
        return  
    return prolog,core,epilog